setwd("E:/5hmc_file/2_5hmc_yjp_bam/ASM/")
file=read.csv("20201120/at.least.one.AShM.in.DC.add.BF.beta0.add.CCHC.csv",head=T)
filea=read.csv("20201112做汇总表/all.FDR.sig.at.least.one.add.direction.same.diff.csv",head=T)
filea$id=paste(filea$Chr,filea$Start,sep = ":")
filea1=filea[filea$FDR.sig>1,]

file$id=paste(file$Chr,file$Start,sep=":")
file1=file[file$pattern.not.rm.dupl.num.DC>1,]
file2=file1[file1$BF_in_DC>1,]
file3=file1[file1$BF_in_DC>10,]
file2=file2[!is.na(file2$unitID),]

sel=grep(names(filea),pattern = "bayes.beta0")
filea$vaf.up.num=rowSums(filea[,names(filea)[sel]]>0,na.rm = T)
filea$vaf.down.num=rowSums(filea[,names(filea)[sel]]<0,na.rm = T)
filea$mean.beta0=rowMeans(filea[,names(filea)[sel]],na.rm = T)
filea$alt.group=ifelse(filea$mean.beta0>0,"up","down")
tmp=data.frame(filea$unitID,filea[,155:158])
write.csv(tmp,"./20210120.H3k.analysis/117012.analysis/117K.ASH.alt.up.or.down.statis.csv",quote=F,row.names = F)

asm=read.table("E:/1.甲基化分析/ASM/ASM_snp-onlyWGS/ASM_log/869727.all.snp.vaf.up.down",head=T,sep="\t")
length(intersect(filea$unitID,asm$rt.unitID))
[1] 44270